N95-11761 


(NASA-CR— 196805) IMPROVED IMAGE 

DECOMPRESSION FOR REDUCED TRANSFORM 
COOING ARTIFACTS (Notre Dame 
Univ.) 10 p 


Unc ! as 


G3/61 0022363 



NASA-CR-1 96805 


/ 

Improved image decompression for reduced 
transform coding artifacts 

Thomas P. O’Rourke, Student Member, IEEE , and 
Robert L. Stevenson, Member, IEEE 



Abstract — The perceived quality of images reconstructed 
from low bit rate compression is severely degraded by the 
appearance of transform coding artifacts. This paper pro- 
poses a method for producing higher quality reconstructed 
images based on a stochastic model for the image data. 
Quantization (scalar or vector) partitions the transform co- 
efficient space and maps all points in a partition cell to a 
representative reconstruction point, usually taken as the 
centroid of the cell. The proposed image estimation tech- 
nique selects the reconstruction point within the quantiza- 
tion partition cell which results in a reconstructed image 
which best fits a non-Gaussian Markov random field (MRF) 
image model. This approach results in a convex constrained 
optimization problem which can be solved iteratively. At 
each iteration, the gradient projection method is used to 
update the estimate based on the image model. In the 
transform domain, the resulting coefficient reconstruction 
points are projected to the particular quantization parti- 
tion cells defined by the compressed image. Experimental 
results will be shown for images compressed using scalar 
quantization of block DCT and using vector quantization of 
subband wavelet transform. The proposed image decom- 
pression provides a reconstructed image with reduced visi- 
bility of transform coding artifacts and superior perceived 
quality. 


I. INTRODUCTION 

Source coding of image data has been a very active area 

— of research for many years. The goal is to reduce the num- 
ber of bits needed to represent an image while making as 
few as possible perceptible changes to the image. Many al- 

gorithms have been developed which can successfully com- 
press a grayscale image to around 0.8 bits per pixel (bpp) 
with almost no perceptible effects [1]. A problem arises, 
however, as these compression techniques are pushed be- 

“ yond this rate. For higher compression ratios (< 0.4 bpp 
for grayscale) most algorithms start to generate artifacts 
which severely degrade the perceived quality of the image. 
The type of artifacts generated is dependent on the com- 
pression technique and on the particular image. For block 
encoded images, the most noticeable artifact is generally 
the discontinuities present at the block boundaries. For 
subband encoded images, it is the ringing at image edges 
which is usually first perceived. 

This paper proposes a technique which addresses this 

— issue through a postprocessing algorithm which greatly re- 
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duces the artifacts introduced by the coding technique. It 
is based on a stochastic framework in which probabilistic 
models are used for both the noise introduced by the coding 
and for a “good” image. The restored image is the maxi- 
mum a posteriori (MAP) estimate based on these models 
[2; 3]. Previous techniques which have tried to address 
this issue have various problems which limit their ability 
to produce high quality image estimates. Some techniques 
propose changes in the way the image is coded[4; 5]; this 
however reduces the efficiency of the source coder and thus 
reduces the compression ratio. The proposed postprocess- 
ing technique only requires modification of the image de- 
coder. Linear based estimators [6], while removing some 
artifacts, usually degrade edge information in the original 
image. Several techniques try to overcome this smoothing 
of the edges by first estimating the edge information in the 
compressed image data [7; 8]. This however is a very diffi- 
cult task for very high compression ratios, where the actual 
edge information is somewhat scrambled. 

This paper will first describe a generic model for im- 
age compression. For the purpose of the reconstruction 
algorithm, this model is descriptive enough to describe 
many compression techniques. A decompression algorithm 
is then described based on a previously proposed image 
model [9; 10]. The computational algorithm is also de- 
scribed. Experimental results are shown for two differ- 
ent transform coding methods. The first method is the 
Joint Photographic Experts Group (JPEG) image com- 
pression standard [11; 1]. The second method uses a sub- 
band wavelet transform with a vector quantization (VQ) 
proposed in [12]. It can be seen that images reconstructed 
with this new postprocessing technique show a reduction 
in many of the most noticeable artifacts and thus allow 
higher compression ratios. 

II. DECOMPRESSION ALGORITHM 

To decompress the compressed image representation, a 
MAP technique is proposed. Let the compressed image 
data be represented by y while the decompressed full reso- 
lution image is represented by z. For MAP estimation, the 
decompressed image estimate z is given by 

z - arg max L(z|y), (1) 

z 

where L(-) is the log likelihood function !(■) = logPr(-) 
which in this case is a measure of how likely the decom- 
pressed image z resulted in the compressed representation 



y. Using Bayes rule 


z — arg max < log 


Pr(y\z)Pr(z) 


z ^ “ Pr( y) ■ ' {2) 

= arg max {log Pr(y |z) + logPr(z) - logPr(y)}(3) 
z 

= arg m|ix {log Pr(y|z) + log Pr(z)} , (4) 


where the Pr(y) term is dropped because it is a constant 
with respect to the optimization parameter z. The condi- 
tional probability Pr(y|z) is based on the image compres- 
sion method while the probability Pr(z) is based on prior 
information about the image data. 

A. Image compression model 

In a transform coding compression technique, a unitary 
transformation H is applied to the original image x. The 
compressed representation y is obtained by applying a quan- 
tization Q to the transform coefficients which can be writ- 
ten as 

y = Q[Hx]. (5) 

Quantization partitions the transform coefficient space and 
maps all points in a partition cell to a representative re- 
construction point, usually taken as the centroid of the 
cell. The indices of these cells are usually entropy coded 
and then transmitted as the compressed representation y. 
In the standard image decompression method, the recon- 
structed image is given by 

z = H~ l Q~ l [y]i (6) 


A Gibbs distribution is used to explicitly write the distri- 
bution of MRF’s. A Gibbs distribution is any distribution 
which can be expressed in the form 



where Z is a normalizing constant, V ,r c (-) is any function of 
a local group of points c and C is the set of all such local 
groups. Note that the Gaussian model is a special case 
of the MRF. To understand how to include discontinuities 
into the statistical model, it is important to first under- 
stand what the model represents and how to define the 
model for a particular application. For a particular source 
signal z\. the value of the probability measure Pr(z\) is 
related to how closely z\ matches our prior information 
about the source. So a Z[ which closely matches our prior 
information should have a higher probability than one that 
does not. For this to be true, the function V r c (-) should pro- 
vide a measure of the consistency of a particular z, where 
a z more consistent with the prior information will have 
smaller values of V c ( )- The situation which is important in 
this work occurs when the prior information is mostly true 
but a limited amount of inconsistency is allowable (e.g. a 
piecewise smooth surface; that is, a surface which is mostly 
smooth but a few discontinuities are allowable). 

In this research a special form of the MRF is used which 
has this very desirable property. This model is character- 
ized by a special form of the Gibbs distribution 


where the inverse quantization maps the indices to the re- 
construction points. 

Since quantization is a many- to-one operation, many im- 
ages map into the same compressed representation. The 
operation of the quantizer is assumed to be noise free; that 
is, a given image z will be compressed to the same com- 
pressed representation y every time. The conditional prob- 
ability for the noise free quantizer can be described by 


Since 


Pr( y|z) 


log Pr(y |z) 


/ l, y = Q[Hz], 
{ 0, y # Q[Hz}. 


/ 0, y = Q[Hz], 

\ -OO, y # Q[Hz], 


(7) 


(8) 


the MAP estimation in (4) can be written as the con- 
strained optimization problem 


z — arg min{— log Pr(z)}, (9) 

where Z is the set of images which compress to y (i.e. 
Z = {z:y = Q[Hz}}). 

B. Image model 

For a model of a “good” image (i.e. Pr( z)) a non- 
Gaussian Markov random field (MRF) model is used [9; 
10]. This model has been shown to successfully model both 
the smooth regions and discontinuities present in images. 


Pr(x) = ^exp{ — j y^pr(dcX)} (11) 

cec 

where A is a scalar constant that is greater than zero, d c 
is a collection of linear operators and the function pr(') is 
given by 

( \ - I u2 ’ M ^ T ’ I'12'l 

| t- +2T(M -t), M>r ; {1 } 

see Figure 1. Since pr() is convex, this particular form of 
the MRF results in a convex optimization problem when 
used in the MAP estimation formulation (4). Therefore, 
such MAP estimates will be unique, stable, and can be 
computed efficiently. The function pr(') is known as the 
Huber minimax function [13] and for that reason this statis- 
tical model is called the Huber-Markov random field model 
(HMRF). 

For this distribution, the linear operators, d c , provide 
the mechanism for incorporating what is considered con- 
sistent most of the time while the function pr(\) is the 
mechanism for allowing some inconsistency. The parame- 
ter T controls the amount of inconsistency allowable. The 
function pr(-) allows some inconsistency by reducing the 
importance of the consistency measure when the value of 
the consistency measure exceeds some threshold, T. 

For the measure of consistency, the fact that the dif- 
ference between a pixel and its local neighbors should be 
small is used; that is, there should be little local variation 
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Fig. 1 Huber minimax function pr( ) superimposed on quadratic func- 
tion. 


in the image. For this assumption, an appropriate set of 
consistency measures is 

{d l c z}c£C — {z m,n - %k j} k , (13) 

where A /* mri consists of the eight nearest neighbors of the 
pixel located at (m, n) and N is the dimension of the image. 
Across discontinuities this measure is large, but the rela- 
tive importance of the measure at such a point is reduced 
because of the the use of the Huber function. 

The MAP estimate can now be written as 


h[i-l] h[i] h[i+I] 

i i i 



Fig. 2. Illustration of projection operator for scalar quantization. 


to find the steepest direction towards the minimum 
g (i) = V (Y,PT(dW ky )) = X>t( d‘z<*>)dj, (17) 

\c£C / c€C 

where p f T (u) is the first derivative of the Huber function. 
The size of the step c*^ is chosen as 

, . , *(*>«*(*) 

a (k) 2 . ( 18) 

“ g (fc)< (Ec€CpT(^ (t) ) d ^)g (fc) ' 

This choice of step size is based on selecting the opti- 
mal step size for a quadratic approximation to the non- 
quadratic function in (15). Since this is an approximation, 
the value of the objective function may increase if the step 
size is too large. To avoid this potential problem, the value 
of a( k) is divided by 2 until the step size is small enough 
that the value of the objective function is decreased. Since 
the updated estimate w (fc+1 \ 

w (t + l) = z r*) + a (4r) g (*) j ( 19 ) 

may fall outside the constraint space 2 , w ( k + l 1 is projected 
onto 2 to give the image estimate at the (fc~ b l)-th iteration 


z = argmmy^V c (z) (14) 

z e2 ' 
c£C 

— = argmin ^ ^ pr(zm y n — zjb.f)- ( 15) 

^ 1 ,n < jV fc ,/£ Af m ,« 

As a result of the choice of image model [9; 10], this results 
in a convex (but not quadratic) constrained optimization 
which can be solved using iterative techniques. 

III. RECONSTRUCTION ALGORITHM 

An iterative approach is used to find z in the constrained 
minimization of (15). An initial estimate is improved 
by successive iterations until the difference between z^ k ^ 
and z^ +1 ^ is below a given threshold e. The rate of con- 
vergence of the iteration is affected by the choice of the ini- 

tial estimate. A better initial estimate will result in faster 

convergence. The initial estimate used here is formed by 
the standard decompression 

z<°> = H~ l Q- l [y\. (16) 

Given the estimate at the Ar-th iteration, z^ k \ the gra- 

— dient descent method is used to find the estimate at the 
next iteration, z^ + 1) . The gradient of ^ prid^z) is used 


z(* +1 > = P*(w<* +l) ). (20) 

In projecting the image w^ +1 ) onto the constraint space 
2 , we are finding the point £ 2 for which ||z^' +1 ^ — 

w^ + 1 )|| is a minimum. If + £ 2, then z £ * +1 ) = 

w^+U and ||z (fc+1 ^ - w^ +1) || = 0. Since H is unitary, 

!|//z ( * +1) - tfw ( * +1) || = ||z (fc+1) - w ( * +1) || (21) 

and the projection can be carried out in the transform do- 
main. 

The form of Vz is dependent on the quantizer Q. The 
projection operators Vz for scalar and vector quantization 
are described in the following subsections. 

A. Scalar quantization 

The scaler quantizer is a partition of the real number line 
7Z and is defined by its breakpoints. The breakpoints are 
determined by the amount of compression desired. The 
boundaries or breakpoints of cell i are l[i ] and h[i], see 
Figure 2. If a coefficient u falls in cell i, 

l[i\ < w < h[i], (22) 

then index i is transmitted. The compressed representation 
y is the set of indices for all the transform coefficients. 
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Projection to the constraint space is rather simple for 
scalar quantization. To find the transform coefficients //. in 
H z( k + l \ find the corresponding coefficient lj in H 
From the compressed representation y, find in which cell 
jjt should fall. Assume /j should fall in cell i. The standard 
decompression method would use (j — r where r is the 
centroid of the cell. The coefficient fi which falls in cell i 
can be determined by 

l[i], u < /[»'], 

w, l[i] <u>< /* [*] , (23) 

h[i], h[i ] < u;. 

This will minimize — # w (*+0|| while assuring 

that y = Q[H The image estimate is ob- 

tained by performing the inverse transformation on H z^' + 1) 

B. Vector quantization 

The vector quantizer is a partition of 7£ n and is defined 
by its codebook vectors r [i] . The size of the codebook is 
determined by the desired bit rate for that subband. 

A coefficient vector uus considered in cell i if it is closer 
to code vector r[z] than to any other codevector 

||u> — t[j]|| < ||u> — r [j ] 1 1 , j*i. (24) 

The distance metric used here is the squared Euclidean 
distance. Although the illustrations in Figure 3 represent 
a 2-D quantizer, the method described below for finding /i 
directly generalizes to higher dimensional spaces. 

Projection to the constraint space is more complicated 
for VQ than for the scalar case. The standard decompres- 
sion method would use [i = r[i] for the reconstruction 
point. To project the vector u; to cell i, it is necessary to 
find the point in cell i which is closest to to The easy case 
occurs when u; itself is in cell i since /x = oj>. In the more 
difficult case where is closer to some other codevector, 
the projection point n will be on the boundary of cell i. 
The first step in projecting u to /x is to move from u to 
the boundary of cell i. The next step is to move along the 
boundary to . 

The points £ on the boundary of cell i all satisfy 


x 



Fig. 3 Illustration of projection operator for vector quantization 


where t 6 [0, 1] if the intersection exists. Among all j ^ i, 
the smallest positive value of t will define the point £ t on 
the cell boundary. £ x will be equally distant from r[i ] and 
r[k}. 

The next step involves traveling along the boundary of 
the cell to /x . The idea is to move towards u> but not leave 
the cell boundary. Let £ 2 be the perpendicular projection 
of onto the k-th hyperplane (||£ — t[z]|| = ||£ - r[£]||) . 
£o is the point in this hyperplane which is closest to and 
can also be found as the intersection of a line 

£ = w 4- t(r[{] - r[k]) (28) 


ll£ ~ r[i]\\ - ||£ - r[j]|| (25) 

for some j ^ i. For a particular value of j, (25) defines 
a hyperplane of points equally distant from the i-th and 
j-th codevectors. Since a; is outside cell i, the line segment 
from u> to t[z] described by 

£ = *u> + (1-*M*1i 0<f<l (26) 

will intersect the boundary of cell i . For each codevector 
r[j ] , j i, it is easy to determine if the intersection of 
the hyperplane (25) with the line segment (26) exists. The 
value of t in (26) which defines the point of intersection is 
given by 

_ ( r H — T [i]) ' ( r [i] — r[J]) 

2(r[i] — w) ■ (r[«] — t[j]) 


with a hyperplane. Both £i and £ 2 will be in this hyper- 
plane. Next, move along the line segment from £ x towards 
£ 2 , stopping at the first intersecting hyperplane. The in- 
tersection point becomes the new £ t and this hyperplane 
becomes the new k - th hyperplane. £ 2 now becomes the 
projection of u; onto this new hyperplane and the process 
is repeated until one of the stopping conditions is met. This 
process is stopped if 

• £j travels all the way to £ 2 without being intercepted 
by another hyperplane 

• or if £ t would leave the cell boundary by moving to- 
wards £ 2 

• or if £ t does not make progress towards to . 

When this process is completed, fi = £ : 
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IV. EXPERIMENTAL RESULTS 

In this section, the results of using the proposed post- 
processing method are shown for images compressed using 
scalar quantization of block DOT and using vector quanti- 
zation of subband wavelet transform. 

A. Scalar quantization 

The standard “Lenna” image shown in Figures 4(a) and 
5(a) was compressed to 0.264 bpp by the JPEG algorithm 
[11; 1]. The result of standard decompression of this image 
is shown in Figure 4(b) and enlarged to show detail in Fig- 
ure 5(b). At this compression ratio, the coding artifacts are 
very noticeable. Most noticeable are the blocking effects of 
the coding algorithm, but aliasing artifacts along large con- 
trast edges are also noticeable. The result of postprocess- 
ing this image is shown in Figure 4(c) and enlarged to show 
detail in Figure 5(c). The number of iterations needed to 
sufficiently reduce the artifacts is dependent on the severity 
of the degradation. More iterations are required for con- 
vergence with more severe degradation. Notice that the 
blocking effects have been removed in the postprocessed 
image. This can be most easily seen in the shoulder and 
background regions. Notice that while the discontinuities 
due to the blocking effects have been smoothed, the sharp 
discontinuities in the original image, such as along the hat 
brim, have been preserved. 

For comparison, the method for reducing blocking arti- 
facts proposed in the JPEG standard [1] is shown in Fig- 
ures 4(d) and 5(d). For a given block, the five lowest 
frequency AC coefficients are predicted from the DC coef- 
ficients of the given block and the 8 nearest blocks. T he 
image is modeled as a 2-D second order polynomial over 
these 9 blocks. This polynomial is defined by requiring that, 
the mean of the polynomial over each block must match 
the DC value for that block. To reduce the blocking, the 
estimated AC coefficients corresponding to a DCT of this 
polynomial are used if they do not contradict the transmit- 
ted information. See [1] for details. Note that while the 
blocking artifacts are reduced by this approach, they are 
still visible. The overall improvement in the reconstructed 
image is much less significant than the improvement seen 
in Figures 4(c) and 5(c). 

While the blocking artifact dominated the standard re- 
construction of the “Lenna^ image, ringing artifacts are 
dominant in document images such as text or graphics [14], 
see Figure 6. While JPEG may not be the most appropri- 
ate compression method for document images, it can be 
seen from the example shown in Figure 6 that ringing arti- 
facts are suppressed by the postprocessing algorithm. The 
postprocessing algorithm removes ringing artifacts as well 
as blocking artifacts. 

B. Vector quantization 

The “Lenna” image was decomposed in the manner of 
Mallat [15] by a three level separable subband structure 
using 2-D separable filters based on 1-D Dan be chics’ 12-th 
order wavelet filters [16]. Each of the subbands is vector 
quantized quantized independently with the transform co- 


efficients grouped into 4 by 4 blocks (or 8 by 8) which are 
treated as 16 (or 64) dimensional vectors. The bit rate 
for each subband is determined by the method given in 
[12] based on the desired overall bit rate and human vi- 
sual system considerations. The bit rate is higher for lower 
frequency subbands. The highest frequency subbands are 
allocated zero bpp and are dropped altogether. For each of 
the remaining subbands, a locally optimal VQ codebook is 
designed for that subband signal. See [17] for general in- 
formation on VQ and [12] for details on the VQ approach 
which is used here. The point here is to demonstrate the 
postprocessing technique is sufficiently flexible to be ap- 
plied to any transform coding method. 

Figure 7 shows the results for moderate compression 
(0.444 bpp) and postprocessing for 9 iterations. The results 
of postprocessing for 20 iterations are shown in Figure 8 for 
the case of heavy compression (0.155 bpp). The greatest 
improvement comes in the first few T iterations. Further it- 
erations of postprocessing continue to improve the image, 
but the amount of improvement in each iteration decreases 
as the image gets better. 

The primary artifacts of this coding technique are alias- 
ing at edges and loss of high frequency information. As ex- 
pected, these artifacts are more severe in the case of heavy 
compression. Most aliasing artifacts disappear after only 
a few iterations in the case of moderate compression. Al- 
though the aliasing artifacts are more persistent for heavy 
compression, the first few iterations still reduce the alias- 
ing a great deal. Even for heavy compression, the aliasing 
artifacts are eventually no longer perceptible. 

Information contained in high frequency subbands, such 
as the details in the hat of the original image, is completely 
lost when those subbands are dropped. The postprocess- 
ing cannot reinvent this information which is completely 
absent from the compressed representation. However, for 
structures such as edges which have information content 
spread across low frequencies as well as high frequencies, 
the postprocessor is able to restore some high frequency 
information based on the low frequency components in the 
compressed representation. This allows crisper edges which 
contain high frequencies to appear in the restored image. 
See for example the hat brim and along the edge of the 
shoulder. As in the scalar quantization case, the struc- 
tures supported by the compressed representation and by 
the image model are restored while artifacts of the trans- 
form coding process are suppressed. 

V. CONCLUSION 

The problem of image decompression has been cast as 
an ill-posed inverse problem, and a stochastic regulariza- 
tion technique has been used to form a well-posed recon- 
struction algorithm. A statistical model for the image was 
produced which incorporated the convex Huber minimax 
function. The use of the Huber minimax function pr{) 
helps to maintain the discontinuities from the original im- 
age which produces high resolution edge boundaries. Since 
Pt( ) is convex, the resulting multidimensional minimiza- 
tion problem is a constrained convex optimization problem. 
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Fig 4. Scalar quantization DOT com| -n-ssr-. I i mag'-s (O.'JiM bpp) (a) original image- 
cessing algorithm, T = 1.0 ( «'l ) I with AO ciH-tficient predict ion 


b) standard decompression (c) After 0-th iteration of postpro- 



Efficient computational algorithms can be used in the min- 
imization. The proposed image decompression algorithm 
produces reconstructed images which greatly reduced the 
noticeable artifacts which exist using standard decompres- 
sion techniques. 
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Fig. 6. Scalar quantization DCT compressed images (0 264 bpp): (a) original image (b) enlarged (c) standard decompression with ringing artifacts (d) 
enlarged (e) After postprocessing algorithm, T = 1.0 (f) enlarged. Note for reviewers: halftoning artifacts distort the fetter edges. This distortion 
is not present in the original images and wilt not appear m the glossy images in the final manuscript. 
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Fig. 7 VQ subband images (0.4-M bpp) (a) standard decompression (b) enlarged (c) after 9-th iteration of postprocessing algorithm, T = 10 (d) 
enlarged. 
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Fig. 8 VQ subband images (0.155 bpp): (a) standard decompression (b) enlarged (c) after 20-th iteration of postprocessing algorithm, T - 1.0 (d) 
enlarged 
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